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ABSTRACT 


In  this  investigation  the  thermal  entrance  region  problem  is  solved 
for  flow  in  rectangular  ducts  of  various  aspect  ratios  and  for  various 
Peclet  numbers.  The  assumptions  under  which  the  problem  is  solved 
are  steady,  fully  developed  laminar  velocity  profile,  constant  fluid  prop¬ 
erties  of  viscosity,  density,  specific  heat  and  thermal  conductivity,  con¬ 
stant  wall  temperature,  and  a  uniform  inlet  fluid  temperature.  Included 
in  the  solution  are  axial  conduction  and  viscous  dissipation.  The  method 
of  B.  G.  Galerkin  is  used  to  formulate  an  approximate  series  solution  of 
the  problem.  The  data  presented  include  bulk  mean  temperature,  local 
Nusselt  number,  and  the  ratio  of  the  local  heat -transfer  rate  to  the  long 
wall  to  the  local  heat -transfer  rate  to  the  short  wall.  It  is  concluded  that 
for  low  Peclet  numbers,  neglecting  the  axial  conduction  term  leads  to  con¬ 
siderable  error  in  the  solution. 


111 


AEDC-TR.69-115 


CONTENTS 

Page 

ABSTRACT .  iii 

NOMENCLATURE .  viii 

I.  INTRODUCTION 

1. 1  Background  Information .  1 

1.  2  Statement  of  the  Problem .  2 

II.  ANALYTIC  PROCEDURE 

2.  1  Formulation  of  the  Problem  with  Governing 

Equations  and  Boundary  Conditions .  2 

2.  2  Nondimensionalized  Equations  .  5 

III.  APPROXIMATE  SOLUTION  OF  THE  ENERGY  EQUATION 

BY  THE  METHOD  OF  B.  G.  GALERKIN .  7 

IV.  NUMERICAL  RESULTS 

4.  1  Numerical  Procedure .  11 

4.  2  Results .  22 

V.  GENERAL  CONCLUSIONS .  58 

REFERENCES .  58 

ILLUSTRATIONS 

Figure 

1.  Flow  Geometry . 3 

2.  Mean  Velocity  versus  Aspect  Ratio .  22 

3.  Bulk  Mean  Temperature  for  Aspect  Ratio  of  10  and 

Parallel  Plates .  23 

4.  Bulk  Mean  Temperature  for  Aspect  Ratio  of  1  ...  .  24 

5.  Bulk  Mean  Temperature  for  Aspect  Ratio  of  1  and 

Various  Peclet  Numbers  . .  25 

6. ‘  Bulk  Mean  Temperature  for  Aspect  Ratio  of  2  and 

Various  Peclet  Numbers .  26 

7.  Bulk  Mean  Temperature  for  Aspect  Ratio  of  5  and 

Various  Peclet  Numbers .  27 

8.  Bulk  Mean  Temperature  for  Aspect  Ratio  of  10  and 

Various  Peclet  Numbers .  '28 

9.  Bulk  Mean  Temperature  for  Pe  of  10  and  Various 

Aspect  Ratios . 29 


v 


AEDC.TR.69.1  15 


Figure  Page 

10.  Bulk  Mean  Temperature  for  Pe  of  100  and  Various 

Aspect  Ratios . 30 

11.  Bulk  Mean  Temperature  for  Pe  of  300  and  Various 

Aspect  Ratios .  31 

12.  Bulk  Mean  Temperature  for  Pe  of  1000  and  Various 

Aspect  Ratios .  32 

13.  Bulk  Mean  Viscous  Temperature  for  Fully  Developed 

Flow  versus  Aspect  Ratio . 33 

14.  Bulk  Mean  Viscous  Temperature  for  Aspect  Ratio  of  1 

and  Various  Peclet  Numbers .  34 

15.  Bulk  Mean  Viscous  Temperature  for  Aspect  Ratio  of  2 

and  Various  Peclet  Numbers .  35 

16.  Bulk  Mean  Viscous  Temperature  for  Aspect  Ratio  of  5 

and  Various  Peclet  Numbers .  36 

17.  Bulk  Mean  Viscous  Temperature  for  Aspect  Ratio  of  10 

and  Various  Peclet  Numbers .  37 

18.  Bulk  Mean  Viscous  Temperature  for  Pe  of  10  and 

Various  Aspect  Ratios .  38 

19.  Bulk  Mean  Viscous  Temperature  for  Pe  of  100  and 

Various  Aspect  Ratios .  39 

20.  Bulk  Mean  Viscous  Temperature  for  Pe  of  300  and 

Various  Aspect  Ratios .  40 

21.  Bulk  Mean  Viscous  Temperature  for  Pe  of  1000  and 

Various  Aspect  Ratios .  41 

22.  Limiting  Nusselt  Number  versus  Aspect  Ratio  for 

Various  Peclet  Numbers  .  . .  42 

23.  Nusselt  Number  for  Aspect  Ratio  of  1  and  Various 

Peclet  Numbers .  43 

24.  Nusselt  Number  for  Aspect  Ratio  of  2  and  Various 

Peclet  Numbers .  44 

25.  Nusselt  Number  for  Aspect  Ratio  of  5  and  Various 

Peclet  Numbers .  45 

26.  Nusselt  Number  for  Aspect  Ratio  of  10  and  Various 

Peclet  Numbers . 46 


AEDC-T  R-69-1 15 


Figure  Page 

27.  Nusselt  Number  versus  Aspect  Ratio  for  Fully 

Developed  Viscous  Temperature .  47 

28.  Viscous  Nusselt  Number  for  Aspect  Ratio  of  1  and 

Various  Peclet  Numbers .  48 

29.  Viscous  Nusselt  Number  for  Aspect  Ratio  of  2  and 

Various  Peclet  Numbers .  49 

30.  Viscous  Nusselt  Number  for  Aspect  Ratio  of  5  and 

Various  Peclet  Numbers .  50 

31.  Viscous  Nusselt  Number  for  Aspect  Ratio  of  10  and 

Various  Peclet  Numbers .  51 

32.  Heat-Transfer  Ratio  Long-to-Short  Wall  for  Aspect 

Ratio  of  2  and  Various  Peclet  Numbers .  52 

33.  Heat-Transfer  Ratio  Long-to-Short  Wall  for  Aspect 

Ratio  of  5  and  Various  Peclet  Numbers .  53 

34.  Heat-Transfer  Ratio  Long-to-Short  Wall  for  Aspect 

Ratio  of  10  and  Various  Peclet  Numbers .  54 

35.  Viscous  Heat-Transfer  Ratio  Long-to-Short  Wall  for 

Aspect  Ratio  of  2  and  Various  Peclet  Numbers  ....  55 

36.  Viscous  Heat-Transfer  Ratio  Long-to-Short  Wall  for 

Aspect  Ratio  of  5  and  Various  Peclet  Numbers  ....  56 

37.  Viscous  Heat-Transfer  Ratio  Long-to-Short  Wall  for 

Aspect  Ratio  of  10  and  Various  Peclet  Numbers  ....  57 


TABLES 


I.  Integer  Pairs  (p,  q)  for  the  Functions  ^i(a,  J3) .  13 

II.  Eigenvalues  and  Integer  Pairs  for  Aspect  Ratio  of  1 

and  Pe  of  1  and  10 .  14 

III.  Eigenvalues  and  Integer  Pairs  for  Aspect  Ratios  of  1 

and  Pe  of  100,  300,  and  1000  .  15 

IV.  Eigenvalues  and  Integer  Pairs  for  Aspect  Ratio  of  1.  5 

and  Pe  of  1000  .  16 

V.  Eigenvalues  and  Integer  Pairs  for  Aspect  Ratio  of  2 

and  Pe  of  10,  100,  300,  and  1000  .  17 


vii 


AEDC.TR.69:il5 


Page 

TABLES  (CONTINUED) 


VI.  Eigenvalues  and  Integer  Pairs  for  Aspect  Ratio  of  5 

and  Pe  of  10,  100,  300,  and  1000  .  18 

VII.  Eigenvalues  and  Integer  Pairs  for  Aspect  Ratio  of  10 

and  Pe  of  10,  100,  300,  and  1000  . .  19 


APPENDIXES 


I. 

SOLUTION  OF  THE  EQUATION  OF  MOTION  FOR  THE 
VELOCITY  PROFILE . 

61 

n. 

SOLUTION  OF  THE  VISCOUS  DISSIPATION  TEMPERA¬ 
TURE  PROFILE . 

64 

HI. 

SOLUTION  OF  THE  INTEGRALS  FOR  Ilf  12,  AND 

13 . 

68 

IV. 

MISCELLANEOUS  INTEGRALS . 

72 

NOMENCLATURE 

Bn 

Eigenvalues  given  in  Eq.  (14) 

CP 

Specific  heat  at  constant  pressure 

h 

Convective  heat -transfer  coefficient 

K 

Thermal  conductivity 

Ll 

Width  of  duct 

l2 

Height  of  duct 

Nu 

Nusselt  number  defined  by  Eq.  (IV-2) 

P 

Pressure 

Pe. 

Peclet  number  defined  by  Eq.  (7) 

Qqt 

Heat-transfer  rate  per  unit  length  to  a  wall 

Qj3 

Heat-transfer  rate  per  unit  length  to  j3  wall 

viii 


AEDC-TR-69.T15 


R  Aspect  ratio  defined  by  Eq.  (7) 

T  Temperature 

T*  Temperature  defined  by  Eq.  (3) 

Tj?  Viscous  dissipation  temperatue 

Tref  Temperature  nondimensionalization  factor  given  by  Eq.  (7) 

Tw  Wall  temperature 

u  Nondimensionalized  fluid  velocity 

um  Nondimensionalized  mean  velocity  defined  by  Eq.  (1-17) 

W  Fluid  velocity 

Wref  Velocity  nondimensionalization  factor  defined  by  Eq.  (7) 

X  Space  coordinate  shown  in  Fig.  1 

Y  Space  coordinate  shown  in  Fig.  1 

Z  Space  coordinate  shown  in  Fig.  1 

a  Nondimensionalized  space  coordinate  defined  by  Eq.  (7) 

(3  .  Nondimensionalized  space  coordinate  defined  by  Eq.  (7) 

6  Nondimensionalized  temperature 

0f  Nondimensionalized  viscous  dissipation  temperature 

0m  Nondimensionalized  bulk  mean  temperature 

A0rp  Nondimensionalized  temperature  difference  between 

inlet  and  duct  wall 

0W  Nondimensionalized  wall  temperature 

Dynamic  viscosity 

?  Nondimensionalized  space  coordinate  defined  by  Eq.  (7) 

p  Density 

<f>i  B.  G.  Galerkin  function  set  defined  by  Eqs.  (26)  and  (27) 

SUBSCRIPTS 

o  Inlet  conditions 

H  Viscous  dissipation  component 

A  Nonviscous  dissipation  component 

00  Fully  developed  conditions 


IX 


A  E  DC-T  R-69-1 15 


SECTION  I 
INTRODUCTION 


1.1  BACKGROUND  INFORMATION 

The  problem  of  the  thermal  entrance  region  heat  transfer  for  flat 
ducts  (parallel  plates)  and  circular  tubes  has  received  considerable 
attention  in  the  literature.  However,  the  similar  problem  for  ducts 
of  finite  aspect  ratio  has  been  investigated  far  less  thoroughly. 

The  heat  transfer  in  laminar  flow  of  an  incompressible  fluid  with ' 
constant  properties  in  the  entrance  region  of  a  circular  tube  was  first 
investigated  analytically  by  Graetz  (Ref.  1)  under  the  assumption  of 
specified  uniform  wall  temperature,  no  axial  conduction,  and  no  viscous 
dissipation.  Prins,  Mulder,  and  Shenk  (Ref.  2)  applied  the  Graetz 
method  to  parallel  plates  under  the  same  assumptions.  Various  authors 
(Refs.  3  through  6)  have  extended  the  problem  to  include  more  complex 
boundary  conditions. 

Sparrow  (Ref.  7)  and  Siegel  and  Sparrow  (Ref.  8)  solved  the  parallel 
plate  entrance  problem  for  simultaneous  development  of  velocity  and 
temperature  profiles  by  use  of  the  von  Karman-Pohlhausen  method  for  a 
range  of  Prandtl  numbers  for  uniform  wall  temperature  and  uniform  heat 
flux,  respectively. 

Schneider  (Ref.  9)  included  axial  conduction  in  the  parallel  plate 
solution  for  slug  flow  of  various  Peclet  numbers  with  finite  wall  resist¬ 
ance  and  both  uniform  and  step  discontinuity  ambient  temperature. 

Hwang  and  Fan  (Ref.  10)  used  finite  difference  techniques  to  solve 
the  problem  of  simultaneous  development  of  velocity  and  temperature 
profiles  for  parallel  plates  for  uniform  wall  temperature  and  uniform 
heat  flux.  Yau  and  Tien  (Ref.  11)  solved  the  same  problem  for  flow  of 
a  non-Newtonian  fluid  with  the  use  of  numerical  techniques. 

Mercer,  Pearce,  and  Hitchcock  (Ref.  12)  presented  experimental 
data  for  laminar  flow  of  air  between  heated  plates  for  Reynolds  numbers 
of  300  to  1500. 

For  the  rectangular  duct,  Clark  and  Kays  (Ref.  13)  used  a  numer¬ 
ical  relaxation  method  to  obtain  limiting  Nusselt  numbers  for  laminar 
flow  in  ducts  with  aspect  ratios  of  1,  2,  and  5  for  constant  heat  flux  and 
aspect  ratios  of  1  and  2  for  constant  wall  temperature  but  without  con¬ 
sideration  of  axial  conduction.  Sparrow  and  Siegel  (Ref.  14)  developed 
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a  variational  method  for  application  to  rectangular  ducts  and  calculated 
the  first  two  eigenvalues  and  eigenfunction  constants  for  the  square  duct 
case  for  constant  heat  flux  but  presented  no  further  information. 

r 

Dennis,  Mercer,  and  Poots  (Ref.  15)  used  a  method  similar  to  the 
Galerkin  method  to  solve  the  constant  temperature  wall  case  for  aspect 
ratios  of  1,  1.  5,  2,  4,  and  8  but  did  not  include  axial  conduction. 

In  this  investigation  the  method  of  B.  G.  Galerkin  as  developed  by 
Eraslan  (Ref.  16)  for  entrance  regions  is  used  to  formulate  approximate 
solutions  for  the  thermal  entrance  region  of  rectangular  ducts  of  various 
aspect  ratios  and  Peclet  numbers.  To  be  included  are  axial  conduction 
and  viscous  dissipation. 


1.2  STATEMENT  OF  THE  PROBLEM 

The  purpose  of  this  investigation  is  to  solve  the  thermal  entrance 
region  problem  for  flow  in  rectangular  ducts  of  various  aspect  ratios 
and  for  various  Peclet  numbers.  The  assumptions  under  which  the  prob¬ 
lem  is  to  be  solved  are  steady,  fully  developed  laminar  velocity  profile; 
constant  fluid  properties  of  viscosity,  density,  specific  heat  and  thermal 
conductivity;  constant  wall  temperatures;  and  a  uniform  inlet  fluid  tem¬ 
perature.  To  be  included  in  the  solution  are  viscous  dissipation  and 
axial  conduction. 

The  assumption  of  fully  developed  laminar  flow  in  a  thermal  entrance 
region  is  applicable  for  fully  developed  flow  in  a  duct  where  the  fluid  en¬ 
counters  a  step  change  in  duct  wall  temperature.  It  may  also  be  applied 
as  an  approximation  for  a  combination  hydrodynamic,  thermal  entrance 
region  where  the  velocity  profile  develops  much  more  rapidly  than  the 
thermal  profile. 


SECTION  II 

ANALYTIC  PROCEDURE 


2.1  FORMULATION  OF  THE  PROBLEM  WITH  GOVERNING  EQUATIONS  AND 
BOUNDARY  CONDITIONS 

The  mathematical  problem  is  formulated  in  the  Cartesian  coordinate 
system  (X,  Y,  and  Z)  with  Z  taken  in  the  direction  of  the  applied  pressure 
gradient,  dP/dZ.  The  origin  of  the  system  is  taken  at  one  corner  of  the 
rectangular  duct  and  in  the  plane  at  which  the  thermal  entrance  region  is 
initiated,  as  shown  in  Fig.  1. 
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For  laminar  incompressible  flow  with  constant  fluid  properties  and 
fully  developed  velocity  profile,  the  equation  of  motion  may  be  written 
as 


with  the  boundary  conditions 

W(D  =  o 


(1) 


where  r  is  the  boundary  of  the  duct.  The  energy  equation  may  be  written 
as 


0  °p w  It  =  K 


32T 


;  ax 


(2) 


with  boundary  conditions 

T(r,Z)  =  Tw 


T(X,Y,0>  =  T  (X,Y) 


Define  a  new  temperature  T*  as 

T*  =  T  -  T. 


=  T(X,Y,Z)  -  Tf(X,Y,Z) 


(3) 


where  Tf  is  the  fully  developed  viscous  temperature  profile.  The  energy 
equation  then  becomes 
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where  T*  satisfies  the  equation 


p  Cn  W 


3T* 

3Z“ 


(Vt*  +  32T*~| 

[ax5”  3Y2  +  3Z2  J 


For  constant  wall  temperature  boundary  conditions, 


3Tf 

3Z“ 


=  0 


and  thus  Tf  satisfies  the  equation 


’32Tf 

32Tf" 

f3W\2  . 

1  3W\2" 

L^~ 

+  — *— 
3Y2 

+  *j[[37')  + 

w| 

The  boundary  conditions  for  T*  and  Tf  are 

t*  (r,z>  =  o 


T*  (X,Y,0)  =  TQ*  (X,Y) 


I  im  T*  =  0 
Z  -*■  » 

Tf  (r,Z)  =  Tw 


(5) 


(6) 


2.2  NONDIMENSIONALIZED  EQUATIONS 


The  mathematical  problem  is  nondimensionalized  based  on  the  width 
of  the  rectangular  duct,  Lj,  and  on  the  pressure  gradient,  dP/dZ,  with 
the  following  variables  and  parameters: 


u 


W 

W 


ref 


e 
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0x  = 


_  'f 


f  T 


ref 


ref 


L  2 

Ll  dP 
u  dZ 


1  ref 


P  wref 
P  cp  Ll 


a  =  X/L, 


6  -  V/Lj 


R  =  L2/Lj 


5  =  Z/L,  Pe 


Pe  = 


P  wref  *-l 


K 


The  variables  or,  j3,  and  ?  are  defined  on  the  intervals 

0  -  a  -  I 

0  -  8  -  R 


0  -  c  -  « 


Substitution  into  Eq.  (5)  for  T*  gives 


u 


80  _  3Z0  320  , 

3a2 


-2  s2e 


— +  Pe 
38  35' 


with  the  boundary  conditions 

0<r,5)  =  0 


0(a,8,O)  =  0Q(a,8) 


(7) 


(8) 


(9) 
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I im  e(a,e,e>  =  0 

£  -*■  00 


Substitution  into  Eq.  (6)  for  Tf  gives 


with  the  boundary  conditions 

ef<r,£)  =  ew 

With  the  substitution 

en  “  ef  "  9' 

Equation  (11)  becomes 


with  the  boundary  condition 

Vr,€>  =  0 


(10) 


(ID 


(12) 


(13) 


SECTION  III 

APPROXIMATE  SOLUTION  OF  THE  ENERGY  EQUATION  BY  THE  METHOD 

OF  B.  G.  GALERKIN 


From  Kantorovich  (Ref.  17)  the  solution  of  an  equation  of  the  form 
L(v)  =  0,  where  L  is  some  linear' differential  operator  in  two  variables 
and  whose  solution  satisfies  homogeneous  boundary  conditions,  may  be 
approximated  by  the  method  of  B.  G.  Galerkin  by  constructing  an  approx¬ 
imate  solution  of  the  equation  in  the  form 


7 


AEDC-TR.69.115 


n 

v(X,Y)  =  ^  a-,  <j>  j  (X,Y) 

1  =  1 


The  (j>  i(X,  Y)(i  =  1,  2, - n)  is  a  chosen  set  of  linearly  independent  func¬ 

tions  satisfying  the  boundary  conditions  and  representing  the  first  n  func¬ 
tions  of  some  system  of  functions  dj^X,  Y)(i  =  1/2, - ®)  which  is  com¬ 

plete  in  the  given  region.  In  order  that  v(X,  Y)  with  n  =  03  be  the  exact 
solution  of  the  given  equation,  it  is  necessary  that  L(v)  be  identically 
equal  to  zero.  This  requirement,  with  L(v)  considered  continuous,  is 
equivalent  to  the  requirement  of  orthogonality  of  L(v)  to  the  functions 
0i(X,  Y).  With  these  conditions  the  system  of  n  linear  equations 


L(v(X,Y) )  $,<X,Y)dX  dY  = 


a  j  <|>j(X,Y)j  (X,Y)dX  dY  =  0 

(I  =  1,2, - n) 


may  be  solved  for  the  n  unknown  coefficients,  Cj,  to  complete  the  approxi¬ 
mate  solution. 


The  assumed  form  of  the  approximate  solution  to  Eq.  (9)  will  be 


K 

0(a,e,O  =  £ 
n=l 


L 


$1 (a,0)e 


"Bn« 


(14) 


where  the  functions  <f>  jfar,  /3)  satisfy  the  boundary  condition 


♦,<r>  =  o 


and  possess  continuous  first-  and  second- order  derivatives  in  a  and  0 
in  the  region  D  and  on  the  closed  boundary,  r,  and  the  a^n^'s  are  con¬ 
stants.  The  constants  a^n)  can  be  determined  by  the  B.  G.  Galerkin 
method. 

Substitution  of  Eq.  (14)  into  Eq.  (9)  and  following  the  method  of 
B.  G.  Galerkin,  which  involves  multiplying  the  equation  by  each  ^^(a,  j3) 
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and  integrating  over  the  region,  D,  yields  the  linear  system  of  equations 


«  0  (15) 


k  =  I ,2,...N 


where  the  integrals  are  given  by 


I  R 


I  I (k, i 


.»■// 


o  o 


V1  <}>,  (a,B)  $k(a,B)  da  d& 


i  « 

l2(k,l)  «  If  u(a,B)  4>.(a,B)  $k(a,B)  da 


dB 


o  o 


R 


1 3< k , T )  b  f  J  $j(a,B>  ♦k(a,B)  da 
Jo  o 


dB 


(16) 


(17) 


(18) 


Since  <b  jfo,  j3)  is  the  known  selected  set,  the  integrals  may  be  eval¬ 
uated  either  analytically  or  numerically  for  a  given  velocity  profile, 
u(or,  |3),  to  give  Eq.  (15)  as  a  system  of  linear  homogeneous  equations 
for  the  undetermined  constants,  a^11'. 

For  a  nontrivial  solution  to  exist  for  the  a-j^'s,  the  determinant 
of  the  coefficients  must  be  identically  zero,  that  is 

Det  0|(k,|>  +  Bn  l2(k,i>  +  (Bn2/Pe2)  l3(k,i)]  =  0  (19) 

Equation  (19)  represents  a  polynomial  of  at  most  2N  degree  in  Bn 
which  must  be  solved  for  the  eigenvalues  <Bn's).  For  each  eigenvalue 
there  will  be  a  set  of  (N-l)  linearly  independent  equations  for  the  ai^n)'s 
which  may  be  solved  in  terms  of  one  of  the  unknown  linear  combination 
constants,  a-j'n*,  which  will  remain  undetermined.  Selecting  the  par¬ 
ticular  undetermined  constant  as  a-^n\  the  linear  system  of  equations 
becomes 
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N 

Y.  c!(n>  I |{k,i)  +  Bn  l2(k,i)  +  (Bn2/Pe2)  i3(k,i)l  = 

1=2  L  J 

-  £l | (k, I )  +  Bn  1 2 C k, I )  +  (Bn2/Pe2)  l3(k, I ) J 

k  =  2,3 ,  • . .N  (20) 


where 


C,(n)  =  a.(n)/a|(n)  (I  *  2) 

With  the  solution  of  the  CjM's  from  the  system  of  Eq.  (20),  the 
approximate  solution  of  Eq.  (9)  becomes 

M  r  N 

=  £  a,(n)  U,(a,B)  +  £  C|(n)  L’8"5  (21) 

n=l  L  i=2  J 

The  constants  a^n)  may  then  be  determined  from  the  inlet  boundary 
conditions  given  for  Eq.  (10), 

M  r  N 

l 

1=2 

which  specifies  the  function  0Q(a,  |8)  as  a  linear  combination  of  known 
functions. 

The  constants  a^n^  may  now  be  determined  by  the  application  of  an 
extension  of  the  Weirstrass  approximation  theorem  (Ref.  20).  Hence, 
multiplication  of  Eq.  (22)  by 


e(a,0,O)=  £  a|(n)  4, ( 

n=  I  L 


a,0)  + 


(n) 


>j(a,6)  = 


0o(a,B) 


(22) 


N 

4, (a, 6)  +  £  Cj(n)  4 j ( a , B ) 

J-2 

and  integration  over  the  region  D  gives 
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M  r  l 

Y  a|(n)  IU,n)  =  JU), 


SL  =  I  ,2,— -M 


(23) 


where  the  integrals  in  Eq.  (23)  are  specified  as 


I  R 


U,n>  =  f  f  a,0)  +  CiCn)  Va'S) 
o  oL  i  =2 


<t> |  (a,S) 


N 

+  Y.  CJ(n) 

J-2 


da  d3 


(24) 


and 


I  R 


o  o 


JU)  J  I  eo(o,0)  *|(a'6)+  Y  Cj(n) 

j"2 


6) 


da  d$ 


(25) 


Since  the  Cj^'s  are  determined  from  Eq.  (20),  Eq.  (23)  represents 
a  nonhomogene ous  system  of  M  linear  equations  in  M  unknowns,  a^11), 
which  has  a  unique  solution  for  j  0o(a,  /3)|  >0.  For  the  thermal  entrance 
region  problem  to  exist,  the  condition  is  satisfied  since  0o(a,  j3)  repre¬ 
sents  the  difference  between  the  duct  inlet  temperature  and  the  fully  de¬ 
veloped  temperature  profiles.  Therefore,  the  a^(n''s  can  be  determined 
to  give  a  complete  approximate  solution  of  the  thermal  entrance  region 
problem. 


SECTION  IV 
NUMERICAL  RESULTS 


4.1  NUMERICAL  PROCEDURE 

The  general  method,  developed  in  the  previous  section  for  Eq.  (9), 
was  used  to  approximate  the  solution  of  the  thermal  entrance  region  heat 
transfer  for  rectangular  ducts  of  aspect  ratio  1,  2,  5,  and  10  and  for 
Peclet  numbers,  based  on  the  pressure  gradient  and  the  velocity  nondi¬ 
mens  ionalization  parameter,  of  10,  100,  300,  and  1000  for  each  aspect 
ratio. 
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The  uniform  wall  temperature  boundary  condition  and  the  resulting 
symmetry  about  the  duct  centerlines  allowed  the  selection  of  the  func¬ 
tions  <b  £)  to  be  of  the  form 

$,(a,8)  ■  sin  p  ira  sin  q  irB/R  (26) 


where  p  and  q  are  pairs  of  odd  integers  corresponding  to  i. 

For  the  special  case  of  the  square  duct,  the  extra  symmetry  of  the 
equal  wall  lengths  allowed  the  functions,  d  i(«,  /3),  to  be  written  in  the 
double  form 


♦  j(o,8)  =  sin  p  if ct  sin  q  ir|5  +  sin  q  wo  sin  p  ttS  (27) 


where  p  and  q  again  are  odd  integers  corresponding  to  i.  This,  in  turn, 
allowed  a  smaller  number  of  functions  to  be  used  for  the  square  duct 
solution  with  a  consequent  decrease  in  computation  time  required.  It 
was  found,  by  solving  specific  cases  with  different  numbers  of  functions, 
that  ten  functions  for  the  square  cases  and  15  functions  for  the  nonsquare 
cases  gave  satisfactory  convergence. 

The  equation  of  motion  (Eq.  (1))  was  solved  in  terms  of  a  double 
Fourier  series  for  the  velocity  profile.  The  solution  is  given  as  Eq. 

(I- 14)  in  Appendix  I.  The  integrals,  Ij(k,  i),  l2(k,i),  and  I3&,  i),  were 
integrated  analytically.  The  integration  is  given  in  Appendix  III.  This 
allowed  a  general  expression  for  the  elements  of  the  determinant  (Eq.  (19)) 
to  be  determined  in  terms  of  the  indices  k  and  i. 

The  determinant  was  found  to  be  diagonally  dominant,  and  the  roots 
of  the  determinant,  Bn's,  were  evaluated  by  use  of  the  maximum  diag¬ 
onal  pivot  Gaussian  reduction  technique  with  the  final  iteration  by  the 
method  of  reguli-falsi.  The  Bn's  were  sufficiently  iterated  to  guarantee 
an  accuracy  of  eight  significant  digits. 

Table  I  gives  the  integer  pairs  (p,  q)  for  the  functions  0i(a,  0).  It 
should  be  noted  that  this  is  a  triangular  truncation  of  the  approximating 
series.  It  was  found  that  this  gave  quicker  convergence  than  a  square 
truncation  for  this  type  of  solution  where  the  coefficients  are  determined 
from  a  set  of  linear  nonhomogeneous  equations.  However,  this  is  not 
necessarily  true  for  a  series  such  as  that  found  for  the  velocity  profile 
solution,  Eq.  (1-14),  in  Appendix  I. 
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TABLE  I 

INTEGER  PAIRS  (p,  q)  FOR  THE  FUNCTIONS  0,  (a,  £) 


1 

R  f  1 

Pe  =  10,  100,  300,  1000 

R  =  1 

Pe  =  100,  300,  1000 

R  =  1 

Pe  -  1 ,  10 

1 

u 

1,1 

1,1 

2 

1.3 

1,3 

1,3 

3 

1,5 

1,5 

1,5 

4 

1,7 

1,7 

1,7 

5 

1,9 

1.9 

1,9 

6 

3,1 

1,11 

1,11 

7 

3,3 

3,3 

3,3 

8 

3,5 

3,5 

3,5 

9 

3,7 

3,7 

3,7 

10 

5,1 

5,5 

3,9 

II 

5,3 

12 

5,5 

13 

7,1 

14 

7,3 

15 

9,1 
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Tables  II  through  VII  give  the  eigenvalues  for  the  cases  of  aspect 
ratios  1,  2,  5,  and  10  and  Peclet  numbers  10,  100,  300,  and  1000.  Also 
included  are  the  cases  of  aspect  ratio  1.  5  and  Peclet  number  1000  and 
aspect  ratio  1  and  Peclet  number  1. 


TABLE  II 

EIGENVALUES  AND  INTEGER  PAIRS  FOR  ASPECT  RATIO 
OF  1  AND  Pe  OF  1  AND  10 


p 

q 

Bn 

Pe  =  1 

Pe  =  10 

1 

i 

4,41418  x  10° 

4. 16395  x  I01 

1 

3 

9.91052  x  10° 

9.69653  x  10* 

3 

3 

1 .33097  x  I01 

1.31405  x  I02 

1 

5 

1.59962  x  I01 

1.57925  x  I02 

3 

5 

1.82987  x  I01 

1.81228  x  I02 

1 

7 

2.21918  x  I01 

2. 19889  x  I02 

3 

7 

2.39072  x  I01 

2.37424  x  I02 

1 

9 

2.84257  x  I01 

2.82231  x  I02 

3 

i 

9 

2,97854  x  I01 

2.96222  x  I02 

1 

1  1 

3.46775  x  I01 

3.44763  x  I02 
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TABLE  III 

EIGENVALUES  AND  INTEGER  PAIRS  FOR  ASPECT  RATIO 
OF  I  AND  Pe  OF  104  304  AND  1000 


p 

q 

Bn 

Pe  =  100 

Pe  =  300 

Pe  =  1000 

1 

i 

2.40433  x  I02 

3.19609  x  I02 

3.36964  x  I02 

1 

3 

7.80820  x  I02 

1.51035  x  I03 

1.94644  x  I03 

3 

3 

1. 15619  x  I03 

2.62382  x  I03 

4.09313  x  I03 

1 

5 

1.38884  x  I03 

3.15950  x  I03 

4.99693  x  I03 

3 

5 

1.64573  x  I03 

4.01269  x  I03 

7.34183  x  I03 

1 

m 

2.00434  x  I03 

4.91449  x  I03 

9.07017  x  I03 

5 

m 

2.05702  x  I03 

5.40188  x  I03 

1.21 179  x  I04 

3 

■ 

2.22453  x  I03 

5.84892  x  I03 

1.32607  x  I04 

1 

2.62970  x  I03 

6.76701  x  I03 

1.40835  x  I04 

1 

1 

3.25869  x  I03 

8.72051  x  I03 

2. 10690  x  I04 
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TABLE  IV 

EIGENVALUES  AND  INTEGER  PAIRS  FOR  ASPECT  RATIO 
OF  l.S  AND  Pe  OF  1000 


p 

a 

Bn 

Pe  '=  1  GOO 

1 

i 

1.76879  x  I02 

1 

3 

7.25395  x  I02 

3 

1 

1.38991  x  I03 

1 

5 

1.79159  x  I03 

3 

3 

2.27422  x  I03 

3 

5 

3.76529  x  I03 

1 

7 

3.82280  x  I03 

5 

1 

4.20636  x  I03 

5 

3 

5.80292  x  I03 

7 

1 

8.21247  x  I03 
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TABLE  V 

EIGENVALUES  AND  INTEGER  PAIRS  FOR  ASPECT  RATIO 
OF  2  AND  Pe  OF  10,  100,  300,  AND  1000 


p 

■ 

u 

Pe  ■ 

175 

Pe  = 

Pe  = 

300 

Pe  = 

TOGO 

1 

1 

3.08423 

X 

I01 

1 .18467 

X 

I02 

1. 31451 

X 

I02 

1.33316 

X 

I02 

1 

3 

5.28847 

X 

I01 

2.96333 

X 

I02 

3.85064 

X 

I02 

4.03407 

X 

I02 

1 

5 

8.09030 

X 

I01 

5.47520 

X 

I02 

8.35539 

X 

I02 

9.24158 

X 

I02 

3 

1 

9. 19682 

X 

I01 

6.54024 

X 

I02 

1 .05383 

X 

I03 

1.19538 

X 

I03 

3 

3 

1.02318 

X 

102 

7.84682 

X 

I02 

1 .38913 

X 

I03 

1 .66049 

X 

I03 

1 

□ 

1.10675 

X 

I02 

8.27962 

X 

I02 

1.43639 

X 

I03 

1 .70509 

X 

I03 

3 

B 

1. 19715 

X 

I02 

9.64982 

X 

I02 

1 .87316 

X 

I03 

2.42366 

X 

I03 

1 

fl 

1. 41147 

X 

I02 

1. 14109 

X 

I03 

2.31374 

X 

I03 

3.21220 

X 

I03 

3 

B 

1. 41921 

X 

I02 

1.201 10 

X 

I03 

2.48366 

X 

I03 

3.27438 

X 

I03 

5 

i 

1 .54327 

X 

I02 

1.25707 

X 

|03 

2.63597 

X 

I03 

4.02052 

X 

I03 

5 

3 

1.60993 

X 

I02 

1 .36135 

X 

I03 

2.87559 

X 

to3 

4. 14845 

X 

I03 

-5 

5 

1 .72741 

X 

I02 

1.51392 

X 

I03 

3.53555 

X 

I03 

6.01251 

X 

I03 

7 

1 

2. 16938 

X 

I02 

1.87356 

X 

I03 

4.14751 

X 

I03 

6.36712 

X 

I03 

7 

3 

2.21957 

X 

I02 

1.98928 

X 

I03 

4.80427 

X 

I03 

8.88796 

X 

I03 

9 

1 

2.79692 

X 

I02 

2.51698 

X 

I03 

6. 16419 

X 

I03 

1 .20650 

X 

I03 
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TABLE  VI 

EIGENVALUES  AND  INTEGER  PAIRS  FOR  ASPECT  RATIO 
OF  5  AND  Pe  OF  10,  100,  300,  AND  1000 


p 

n 

Pe  =  10 

Pe  =  100 

Pe  =  300 

Pe  =  1000 

1 

i 

2.71323  x  I01 

8.82645  x  I01 

9.45499  x  I01 

9.53839  x  IO1 

1 

3 

3.19397  x  I01 

1 . 18039  x  I02 

1.29705  x  I02 

1 .31344  x  IO2 

1 

5 

3.98096  x  I01 

1.70549  x  I02 

1.95330  x  I02 

1.99152  x  IO2 

1 

B 

4.94473  x  I01 

2.41554  x  I02 

2.91486  x  I02 

3.00176  x  IO2 

1 

KB 

6.01 129  x  I01 

3.35243  x  I02 

4.4131 1  x  I02 

4.64455  x  IO2 

3 

9.031 14  x  lO1 

6.10023  x  |02 

9.30171  x  I02 

1.02871  x  IO3 

3 

n 

9.21879  x  lO1 

6.34871  x  I02 

9.84766  x  I02 

1.09692  x  IO3 

3 

9.55470  x  lO1 

6.75797  x  I02 

1.08452  x  IO3 

1.22914  x  IO3 

3 

1.00300  x  I02 

7.47551  x  I02 

1.31498  x  I03 

1.57715  x  IO3 

5 

i 

1.53107  x  I02 

1.21204  x  I03 

2.29434  x  I03 

2.90775  x  IO3 

5 

3 

1.54326  x  I02 

1.23440  x  I03 

2.38673  x  I03 

3.09339  x  IO3 

5 

5 

1.56509  x  I02 

1.29307  x  I03 

2.67877  x  I03 

3.77781  x  IO3 

7 

1 

2.15909  x  lO2 

1.82828  x  I03 

3.90311  x  I03 

5.66938  x  IO3 

7 

3 

2.16910  x  I02 

1.8731 1  x  lO3 

4.24036  x  iO3 

6.97333  x  IO3 

9 

1 

2.78759  x  I02 

2.46842  x  I03 

5.89038  x  IO3 

1.09241  x  IO4 
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TABLE  VII 

EIGENVALUES  AND  INTEGER  PAIRS  FOR  ASPECT  RATIO 
OF  10  AND  Pe  OF  10,  100,  300,  AND  1000 


p 

■ 

Bn 

TO 

Pe  = 

100 

Pe  = 

300 

Tooo 

1 

1 

2.66084 

X 

I01 

8.49501 

X 

I01 

9.07163 

X 

I01 

9.14761 

X 

I01 

1 

3 

2.78789 

X 

I01 

9.20748 

X 

I01 

9.88936 

X 

I01 

9.98038 

X 

I01 

1 

5 

3.02449 

X 

I01 

1 .05706 

X 

I02 

1 .14767 

X 

I02 

1.16006 

X 

I02 

1 

B 

3.3491 | 

X 

I01 

1.25430 

X 

I02 

1.38313 

X 

I02 

1 .40138 

X 

I02 

1 

B 

3.74307 

X 

I01 

1.52746 

X 

I02 

1.72788 

X 

I02 

1.75807 

X 

I02 

3 

B 

9.01022 

X 

I01 

6.06540 

X 

I02 

9.22175 

X 

I02 

1 .01866 

X 

I03 

3 

B 

9.05687 

X 

I01 

6.11968 

X 

I02 

9.33374 

X 

I02 

1 .03235 

X 

I03 

3 

B 

9. 14621 

X 

I01 

6.23653 

X 

I02 

9.62644 

X 

I02 

1.07083 

X 

I03 

3 

9.27769 

X 

I01 

6.51246 

X 

I02 

1 .04882 

X 

I03 

1 . 19224 

X 

I03 

5 

i 

1 .52964 

X 

I02 

1.20916 

X 

I03 

2.28500 

X 

I03 

2.89177 

X 

I03 

5 

3 

1 . 53264 

X 

I02 

1.21497 

X 

I03 

2.31009 

X 

I03 

2.94085 

X 

I03 

5 

5 

1 .53871 

X 

I02 

1 . 23987 

X 

I03 

2.46839 

X 

I03 

3.33297 

X 

I03 

7 

1 

2.15795 

X 

I02 

1 .82484 

X 

I03 

3.88155 

X 

I03 

5.60739 

X 

I03 

7 

3 

2. 16058 

X 

I02 

1 .84707 

X 

I03 

4.1 171  1 

X 

I03 

6.59961 

X 

I03 

9 

1 

2.78647 

X 

I02 

2.46319 

X 

I03 

5.86196 

X 

I03 

1 .08121 

X 

I04 
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For  the  lower  Peclet  numbers  the  integrals  I^(k,  i)  and  Ig(k,  i),  as 
specified  by  Eqs.  (16)  and  (18),  respectively,  dominate  the  determinant 
to  the  extent  that  two  integer  pairs  with  the  same  value  of  the  sum 
(p2  +  q2/R2)  produced  nearly  equal  eigenvalues.  For  the  cases  of  aspect 
ratio  1  and  Peclet  numbers  1  and  10,  the  eigenvalues  produced  by  the 
pairs  (1,  7)  and  (5,  5)  differed  only  in  the  fifth  significant  digit.  This  was 
close  enough  to  produce  difficulty  in  evaluating  the  nonhomogeneous  sys¬ 
tem,  Eq.  (23),  and  as  shown  in  Table  I  for  these  cases,  the  integer  pair 
(3,  9)  was  substituted  for  the  pair  (5,  5). 

It  was  found  that  each  function  i(a,  /3)  produced  an  eigenvalue,  and 
for  the  lower  Peclet  numbers  it  was  possible  to  identify  which  function 
generated  which  eigenvalue  on  the  basis  of  the  value  of  p2  +  q2/R2.  In 
Tables  II  through  VII  the  eigenvalues  are  listed  in  ascending  order. 
Listed  with  the  eigenvalues  are  the  integer  pairs  which  generated  them. 

The  inlet  temperature  boundary  condition  90(a,  f$),  which  was  defined 
as  the  difference  between  the  temperature  at  the  duct  inlet  and  the  final 
developed  viscous  temperature  profile,  was  written  as 

60<a,B)  =  A0t  -  6n(a,B)  (28) 


where  A 0T  is  the  temperature  difference  between  the  constant  tempera¬ 
ture  inlet  condition  and  the  uniform  wall  temperature,  and  0^(tt,  j3)  is  the 
solution  of  the  transformed  viscous  dissipation  equation  given  as  Eq. 

(II- 10)  in  Appendix  n. 

The  constants  a^n^  from  Eq.  (23)  were  solved  for  by  the  maximum 
diagonal  pivot  Gaussian  reduction  technique.  With  separation  of  the  in¬ 
let  boundary  conditions  into  nonviscous  and  viscous  dissipation  compon¬ 
ents,  two  sets  of  constants,  a^x^  and  alr^n^»  were  solved  for  during 
the  same  matrix  reduction.  This  was  accomplished  by  writing  Eq.  (25) 
for  J(i)  as 

JU)  =  JAU>  -  Jn(£>  (29) 


where 


(a, 8) 


+ 


N 


I 


■ 


da  d8 


(30) 
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Y*> 


1  r  r  n  1 

-/  J  0^0,6)  j  |  ( a , 3 )  +  Cj*n* 

0  0  n  L  j=2 


da  dg 


■//  1  1  Ts(s,+)  sin  sira  sin  t  irg/R 

0  0  s=l ,3,5. .  t= I , 3 , 5 . . 


N 

<> I  Ca , P )  +  ^  Cj*n* 

j=2 


and  by  specifying  J^U)  and  J^fi)  as  separate  column  vectors  in  the 
matrix  solution  for  ai^n)  and  Given  in  Appendix  IV  as  Eqs. 

(IV- 1)  through  (IV- 9)  are  the  integrated  solution  forms  for  Eq.  (24) 
for  I(j2,n)  and  Eqs.  (30)  and  (31)  for  JxU)  and  J^U). 

The  temperature  solutions  now  become,  for  the  nonviscous  com¬ 
ponent, 

M 

0x(a,6,5)  =  ASy  aIX(0> 

n=i 

N 

+  Z  C,(n)  <M(a,B) 

1=2 

and  for  the  viscous  component, 

M 

en(a,8,5>  =  Pe  Y*  “It/"1 

n=  I 


4>  |  (a,  B)  + 
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N 


1=2 


4 j (a, 6) 


“Bp* 


(33) 


4.2  RESULTS 

Of  primary  interest  in  this  study  are  the  bulk  mean  temperature, 
0m,  as  specified  by  Eqs.  (IV- 10)  through  (IV- 14)  in  Appendix  IV,  the 
Nusselt  numbers,  Nu,  as  specified  by  Eqs.  (IV-20)  and  (IV-21)  in 
Appendix  IV,  and  the  ratio  of  the  long-to- short  wall  heat-transfer  rates, 
Qj3 lQa>  as  may  be  determined  from  Eqs.  (IV- 15)  through  (IV-19)  in 
Appendix  IV. 

Figure  2  gives  the  nondimens ionalized  mean  velocity,  um,  as  eval¬ 
uated  from  Eq.  (1-17),  as  a  function  of  aspect  ratio.  It  is  this  mean  ve¬ 
locity  which  must  be  multiplied  by  the  Peclet  number  used  in  this  study 
to  obtain  the  more  standard  Peclet  number. 


Fig.  2  Mean  Velocity  versus  Aspect  Ratio 
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Figure  3  gives  a  comparison  of  the  nonviscous  component  of  the 
bulk  mean  temperature  for  aspect  ratio  of  10  and  Peclet  numbers  of 
1000  and  100  to  that  given  by  Ref.  2  for  the  parallel  plate  solution  with 
no  axial  condition.  The  parameter,  Pe',  is  the  Peclet  number  based  on 
the  mean  velocity  for  parallel  plate  flow  as  used  by  Ref.  2.  Figure  4 
gives  a  comparison  of  the  nonviscous  component  of  the  bulk  mean  tem¬ 
perature  for  aspect  ratio  of  1  and  Peclet  numbers  of  1000  and  100  to 
that  given  by  Ref.  15  for  no  axial  conduction.  For  aspect  ratio  of  1, 
the  Peclet  number  of  1000  as  used  in  this  study  corresponds  to  a 
Peclet  number,  based  on  the  mean  velocity,  of  35-.  The  data  from 
Ref.  15  should  be  applicable  for  Peclet  numbers,  based  on  the  mean 
velocity,  of  100  and  above.  It  may  be  seen  from  this  figure  that  ignor¬ 
ing  axial  conduction  for  low  Peclet  numbers  can  give  considerable  error. 


Fig.  3  Bulk  Mean  Temperature  for  Aspect  Ratio  of  10  and  Parallel  Plates 
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Fig.  4  Bulk  Mean  Temperature  for  Aspect  Ratio  of  1 
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Figures  5  through  8  give  the  nonviscous  component  of  the  bulk  mean 
temperature  for  aspect  ratios  of  1,  2,  5,  and  10,  respectively,  with  the 
Peclet  number  as  a  parameter.  Figures  9  through  12  give  the  same  tem¬ 
perature  distribution  for  Peclet  numbers  10,  100,  300,  and  1000,  respec¬ 
tively,  with  the  aspect  ratio  as  a  parameter. 


0  1.0  2.0  3.0  4.0  5.0 


Z/L, 

Fig.  5  Bulk  Mean  Temperature  for  Aspect  Ratio  of  1  and  Various  Peclet  Numbers 
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Fig.  7  Bulk  Mean  Temperature  for  Aspect  Ratio  of  5  ond  Various  Peclet  Numbers 
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Fig.  10  Bulk  Mean  Temperature  for  Pe  of  100  and  Various  Aspect  Ratios 
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Fig.  11  Bulk  Mean  Temperature  for  Pe  of  300  and  Various  Aspect  Ratios 
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Fig.  12  Bulk  Mean  Temperature  for  Pe  of  1000  and  Various  Aspect  Ratios 
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Figure  13  gives  the  bulk  viscous  dissipation  mean  temperature  for 
fully  developed  flow  as  a  function  of  aspect  ratio.  Included  is  the  value 
for  the  parallel  plate  solution. 


Fig.  13  Bulk  Mean  Viscous  Temperature  for  Fully  Developed  Flow  versus  Aspect  Ratio 
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Figures  14  through  17  give  the  viscous  component  of  the  bulk  mean 
temperature  for  aspect  ratios  of  1,  2,  5,  and  10,  respectively,  with  the 
Peclet  number  as  a  parameter.  Figures  18  through  21  give  the  same 
temperature  distribution  for  Peclet  numbers  of  10,  100,  300,  and  1000, 
respectively,  with  the  aspect  ratio  as  a  parameter. 
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Fig.  16  Bulk  Mean  Viscous  Temperature  for  Aspect  Ratio  of  5  and  Various  Peclet  Numbers 
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Fig.  18  Bulk  Mean  Viscous  Temperature  for  Pe  of  10  and  Various  Aspect  Ratios 


38 


lO3/^ 


AEDC-TR-69-1 15 


5.0 


4.0 


3.0 


x 

fir 

<de 

2.0 


1.0 


Fig.  19  Bulk  Mean  Viscous  Temperature  for  Pe  of  100  and  Various  Aspect  Ratios 
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Fig.  20  Bulk  Mean  Viscous  Temperature  for  Pe  of  300  and  Various/Aspect  Ratios 
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Fig.  21  Bulk  Mean  Viscous  Temperature  for  Pe  of  1000  ond  Various  Aspect  Ratios 
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Figure  22  gives  the  limiting  Nusselt  number  as  a  function  of  aspect 
ratio  for  the  various  Peclet  numbers,  which  is  .identical  for  both  the 
viscous  and  nonviscous  components.  Plotted  in  the  figure  are  data  from 
Refs.  13  and  15.  The  values  from  Ref.  13  are  2.  89  and  2.  54  for  the 
aspect  ratios  of  1  and  2,  respectively,  and  from  Ref.  15  are  2.  98  and 
2.  54.  These  values  compare  well  with  the  values  obtained  in  this  study 
for  Peclet  number  of  1000,  which  are  2.  981  and  2.  548.  , 


R 


Fig.  22  -Limiting  Nusselt  Number  versus  Aspect  Ratio  for  Various  Peclet  Numbers 
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give  the  Nusselt  number  for  the  nonviscous 
distance  from  the  duct  entrance  for  aspect 
respectively,  with  the  Peclet  number  as  a 


Pe  *  1000 
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Fig.  24  Nusselt  Number  for  Aspect  Ratio  of  2  and  Various  Peclet  Numbers 
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Fig.  25  Nusselt  Number  for  Aspect  Ratio  of  5  and  Various  Peclet  Numbers 
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Fig.  26  Nusselt  Number  for  Aspect  Ratio  of  10  and  Various  Pec  let  Numbers 
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Eraslan  and  Snyder  (Ref.  18)  have  shown  that  for  viscous  dissipa¬ 
tion  the  mean  value  of  the  dissipation  function  over  the  duct  cross  sec¬ 
tion  is  equal  to  the  product  of  the  pressure  gradient  and  the  mean  ve¬ 
locity.  With  the  values  of  bulk  mean  temperature,  mean  velocity,  and 
Nusselt  number  for  fully  developed  viscous  flow  between  parallel  plates 
obtainable  from  a  closed  form  solution,  it  is  possible  to  obtain  a  more 
accurate  estimation  of  the  fully  developed  viscous  dissipation  Nusselt 
number  for  various  aspect  ratios  by  means  of  the  ratio  of  energy  dissi¬ 
pated  than  is  obtainable  by  differentiating  the  finite  series  for  the  fully 
developed  viscous  temperature  profile.  The  Nusselt  number  may  then 
be  shown  to  be 


% 


Nu 


um  9mnf 
emn  umf 


(34) 


where  the  subscript  f  denotes  parallel  plate  values  and  where  um  and 
0mr )  may  be  obtained  from  Figs.  2  and  13.  Figure  27  gives  the  Nusselt 
number  obtained  by  this  approach  as  a  function  of  aspect  ratio. 


o  2  4  6  8  10 


Fig.  27  Nusselt  Number  versus  Aspect  Ratio  for  Fully  Developed  Viscous  Temperature 
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Figures  28  through  31  give  the  Nusselt  number  for  the  viscous  com¬ 
ponent  as  a  function  of  distance  from  the  duct  entrance  for  aspect  ratios 
of  1,  2,  5,  and  10,  respectively,  with  the  Peclet  number  as  a  parameter. 
These  figures  show  that  the  Nusselt  number  for  the  viscous  component 
is  not  infinite  at  the  duct  entrance,  as  in  the  case  for  the  nonviscous 
component,  but  is  given  by  Fig.  27.  Therefore,  the  curves  in  Figs.  28 
through  31  were  started  at  the  duct  entrance  using  the  values  given  by 
Fig.  27. 


Fig.  28  Viscous  Nusselt  Number  for  Aspect  Rotio  of  1  and  Various  Peclet  Numbers 
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Fig.  29  Viicous  Nusselt  Number  for  Aspect  Ratio  of  2  and  Various  Peclet  Numbers 
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Fig.  30  Yiscous  Nusselt  Number  for  Aspect  Ratio  of  5  and  Various  Peclet  Numbers 
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Fig.  31  Viscous  Nusselt  Number  for  Apsect  Ratio  of  10  and  Various  Peclet  Numbers 
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Figures  32  through  34  give  the  ratio  of  heat-transfer  rate  long-to- 
short  wall  as  a  function  of  distance  from  the  duct  entrance  for  the 
nonviscous  component  and  for  aspect  ratios  of  2,  5,  and  10,  respectively, 
for  the  various  Peclet  numbers.  Near  the  entrance  to  the  duct  the  solu¬ 
tion  converges  more  rapidly  to  the  long  than  to  the  short  wall.  For  this 
region  the  heat-transfer  rate  ratio  has  been  estimated  by  extending  the 
curves  as  dashed  lines. 


Fig.  32  Heat-Transfer  Ratio  Long-to-Short  Wall  for  Aspect  Ratio  of  2  and  Various  Peclet  Numbers 
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Figures  35  through  37  give  the  ratio  of  heat-transfer  rate  long-to- 
short  wall  as  a  function  of  distance  from  the  duct  entrance  for  the  vis¬ 
cous  component  and  for  aspect  ratios  of  2,  5,  and  10,  respectively,  for 
the  various  Peclet  numbers.  Again,  near  the  entrance  of  the  duct,  the 
ratios  have  been  estimated  by  dashed  lines. 


Z/L, 


Fig.  35  Viscous  Heat-Transfer  Ratio  Long-to-Short  Wall  for  Aspect  Ratio  of  2  and  Various  Peclet 
Numbers 
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Fig.  36  Viscous  Heot- Transfer  Ratio  Long-to-Short  Wall  for  Aspect  Ratio  of  5  and 
Various  Peclet  Numbers 
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Viscous  Heat>Transfer  Ratio  Long-to-Short  Wall  for  Aspect  Ratio  of  10  and  Various 
Pec  let  Numbers 
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SECTION  V 

GENERAL  CONCLUSIONS 


It  is  concluded  that  the  B.  G.  Galerkin  method  may  be  used  to  solve 
the  thermal  entrance  region  heat-transfer  problem  for  laminar  flow  of  a 
constant  property  fluid  in  rectangular  ducts  with  axial  conduction  and 
viscous  dissipation  included  if  consideration  is  given  to  avoiding  near 
multiple  roots  for  the  eigenvalues. 

The  separation  of  the  inlet  boundary  conditions  into  viscous  dissipa¬ 
tion  and  nonviscous  dissipation  components  allowed  a  solution  to  be  ob¬ 
tained  which  was  general  for  the  duct  entrance  to  wall  temperature  dif¬ 
ference,  p.  The  combination  of  the  viscous  and  nonviscous  compon¬ 
ents  gives  the  solution  for  the  inlet  boundary  conditions  as  stated  in 
Section  1.2.  The  nonviscous  component  alone  will  give  the  solution  for 
no  viscous  dissipation  and  uniform  temperature  inlet,  as  well  as  for  fully 
developed  viscous  temperature  flow  with  a  step  change  in  wall  tempera¬ 
ture  in  the  Z-direction.  The  viscous  component  alone  will  give  the  solu¬ 
tion  of  viscous  temperature  flow  for  the  condition  of  uniform  inlet  tem¬ 
perature  equal  to  a  uniform  wall | temperature. 

The  solutions  for  the  thermal  entrance  region  in  rectangular  ducts 
with  constant  wall  temperature  conditions  verify  the  fact  that  for  low 
Peclet  numbers,  the  few  previous  cases  solved  for  no  axial  conduction 
can  be  in  considerable.  error. 
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APPENDIX  I 

SOLUTION  OF  THE  EQUATION  OF  MOTION  FOR  THE  VELOCITY  PROFILE 


Given  the  equation  of  motion  (Eq.  (1)) 


a2W  a2W  \_  dP 
v  d? 


ax'-  aY' 


(i-i)' 


and  nondimensionalizing  with 


ref 


Ll2  dP  _  W 
v  dZ  “  u 


(1-2) 


and 


a  & 


(1-3) 


the  equation  becomes 


Wref 

17“ 


a2u  a2u  ' 


I  dp 

p  dZ 


(1-4) 


or 


a2u  a2p 

7?  *7? 


=  -  i 


(1-5) 


The  boundary  conditions  are 


u(D  =  0 


(1-6) 


The  solution  may  be  obtained  using  the  finite  Fourier  sine  transform 
as  given  by  Sneddon  (Ref.  19).  Taking  the  transform  of  both  sides  of  the 
equation 
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i  r\ 

II 

o  o 


Vu  a^u 

1 3a2  332 


r  •  niT®  j 

sin  m  Iras  in  — p  da  d0  = 


I  K 

-// 

o  o 


sin  m  ira  sin  JIIL®.  da  dB 
R 


(1-7) 


the  transform  of  the  left  side  of  the  equation  with  the  specified  boundary- 
conditions  becomes 


-IT 


("?  +  i?)  T  Mm'n> 


(1-8) 


The  transform  of  the  right  side  becomes 

R 


mnit 


£cos  mir  -  ij  cos  mr  -  I 


(1-9)  i 


R 


mnir 


C2D  121  for  m  and  n  odd 


=  0  for  m  or  n  even 


(I- 10) 

(I- 11) 


Then 


Atm 


,  i6  r 

»n)  =  T 


^  R^m^n+mn^  J 


(1-12) 


The  solution  then  becomes 


E  E 


u  =  L-t  L-.  A(m,n)  sin  m  ira  sin  Uli  (1-13) 

m=l,3,5...  n=l,3,5... 


or 
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“  -  Z  Y.  "4  fiAW 

m=l,3,5...  n=l,3,5...  L  ■* 


sin  m  ira  sin  ^5—  (1-14) 


The  mean  velocity  may  then  be  determined  as 

1  R 


-■Jh 


dad  (J/R 


o  o 


(1-15) 


or 


I  R 


If16 

R2 

sin  m  Trn  sin  n7r®  dnflP 

J  J  „4 
o  o  * 

R^m^n+mn^ 

R 

Uli  dodB  (1-16) 


m 


R 


Then 


E  Z 

m- I ,3,5. • .  n™l,3,5«.« 


R^m^n^+m^n^ 


(1-17) 
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APPENDIX  II 

SOLUTION  OF  THE  VISCOUS  DISSIPATION  TEMPERATURE  PROFILE 


Equation  (13)  for  the  viscous  dissipation  temperature  profile 


^  ♦!% 
3  a  38^ 


(II- 1) 


with  boundary  conditions 


=  0 


may  be  solved  using  the  finite  Fourier  sine  transform.  Taking  the 
transform  of  both  sides  as 


R 


II 

o  o 


3"  P  +  I  sin  p  ira  sin  q  dad8 

382  /  R 


3a 


2 


-  Pe 


}/' 

o  o 


3u 

3a 


3u 

38 


sin  p  ira  s*n  q  JliL  dadS 
R 


(II- 2) 


the  transform  of  the  left  side  as  given  by  Sneddon  (Ref.  19)  becomes 

,2 


-  Tt 


2  .  q‘ 


P  + 


rZ 


(H-3) 


From  the  velocity  solution 


3 u 
3a 


2  _  \  Y  V  V  Y 

i - 

CM 

cn 

i _ 

ttRZ 

'M  l  L  L  L 

1  1  m=l ,3,5  n= 1 ,3,5  s=l,3,5  t=l,3,5 

R2m^n+n^J 

R2s2t+t^ 

cos  m  ira  sin  n  !l  cos  s  ira  sin  t 

R  R 


(II- 4) 
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3u\2  1 6 

30  +  -T 


l  E  I  I  I: 

=1,3,5  n= 1,3,5  s=l,3,5  t=l,3,5L 


ttR  ttR 

R2m3+mn2  R2s^+st2 


sin  m  tto  cos  n  lA  sin  s  ira  cos  t  ijL 
R  R 


The  right  side  then  becomes 


„  I6\2 

-Pe  — 


HZ 


m  n  s 


TT2  R4 


I  (R2m2n+n3)(R2s2t+t3) 


(n-5) 


I  K 

I  /  cos  m  ira  cos  s  ira  sin  p  ira  sin  n  2lB  s  i  n  t  s  i  n  q 
J  J  R  R  ^ 


dad  8  + 


III  If _ ££ _ 

m  n  s  t  (R2m5+mn2)(R2s5+s+2) 


i  n 

II 


sin  m  ira  sin  s  na  sin  p  Tra  cos  n  IJL  cos  t  JlE  sin  q  JL§.  dad8 

R  R  R 


j  (H-6) 


By  use  of  trigonometric  identities  for  multiple  angles  the  first  inte¬ 
gral  becomes 


L  R  r 

U\* 

o  o  L 


sin  (m  +  (s+p))ira  -  sin  (m  -  (s+p))ira 


-  sin  (m  +  (s-p))ira  +  sin  (m  -  ( s-p ) ) ttoJ  |^s i n  (n  +  (t-q))!|- 

+  sin  (n-(t-q) )lJL  -  sin  (n+(++q)  )Hjl  -  sin  (n-(t+q) )£§. 

R  R  R 


da  d8 


(II- 7) 


65 


AEDC-TR.69-n5 


which,  on  integration,  becomes 


R  _ I 

rr  m+(s- 


I _ I  +  _ 

m-Ts+pT  m+(s-p)  m 


-  ' I 

-<s-p)J 


(n-8) 


■  ,  J  i  i 

m+(t-q)  n-(t-q)  n+(t+q)  ”  n-(t+q! 


and  exists  only  for  p  and  q  =  1,3,5 - since  m,  n,  s,  and  t  are  odd 

and 


sin  n  irotda  =  0  for  n  even 


By  the  same  process  the  second  integral  becomes 


R  I  I _ I 

4ir2  I  m+(s-p)  m-(s-p)  m+(s 


n+(t+q)  "  ri-(++q)  n+(t-q) 


I  1 

n-(t-q) J 


(II- 9) 


and  exists  only  for  p  and  q  =  1,  3,  5 


The  final  solution  then  becomes 


6n(a,8)  =  Pe  ^  Z-  fs(p,q> 

p-1,3,5..  q-  1 ,3,5. . 


*  fe^P»< 


stn  p  ua  sin  q 


(II- 10) 


where 


-  _  256  R4 

fs<p'q)  ~  “  10 


R2p2+q2 


iii? 
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APPENDIX  111 

SOLUTION  OF  THE  INTEGRALS  FOR  h,  l2,  AND  l3 


For  approximating  functions  of  the  form 

=  sin  p  ira  sin  q  ZJL 
R 


where  p  and  q  are  functions  of  i  and  for  the  velocity  solution 


L.  Z. 


4  m=l ,3,5. . .  n=l,3>5...\R2m3n+mn3 


sin  m  Ttosin  n 


itB 


and 


(k,i)  =  J  J 


-  j  j  $.(a,B)  ^(a^S)  dadS 

o  o 


V2  =  -  *2  (p2  +  q2/R2)  sin  p  ira  sin  q 


then 


Ik 


i  rs 

.»■// 


-  ir2  (p2  +  qz/R2)  sin  p  ira  sin  s  ira  sin  q 


o  o 


sin  t  HjL  dadS 


where  p  =  s  and  q  =  t  for  k  =  i 


l|(k,i)  =  -  n2(p2  +  q2/R2)  £  for  k  = 

4 


I  j  (k,  I )  =  0  for  k  t  i 


(III-l) 

(HI- 2) 

(III- 3) 

(III-4) 

(III- 5) 

(III- 6) 
(III- 7) 
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i  r\ 

1 2<  k,  i )  =  J  J  u<a,0)  <f>j(a,0)  <}>k<a,0>  da  d0 


(III- 8) 


then 


1 2 f  k,  i ) 


■3  II  I 


o  o  m=|,3,5  n=  1,3,5  l^m^n+mn-5 


sin  m  na  sin  p  na  sin  s  ira  sin  nJL®.  sin  sin  +  2L^.da  dB  (III- 9) 

R  R  R 


which,  with  the  use  of  trigonometric  identities,  becomes 


l2(k,i)  «  ^  /  /  Y.  E  16  rVWv5) 

o  o  rn—  i  9 3 9 5 •  •  n_  i  ^3^5*  •  i  * 

sin  (m+(p-s))na  +  si n(m-<p-s) ) na  -  si n(n+f p+s) )na  -  si n(m-(p+s) )na 

sin  (n+(q-t) )  *£.  +  sin  <n-(q-t)l  ij  -  sin  (n+(q++)^. 

R  R  R 


sin  (n-(q++))  -7^  da  dB 
K 


(in- 10) 


Since  p,  q,  s,  and  t  may  be  specified  as  odd  integers  because  of 
symmetry  about  the  duct  centerlines,  the  solution  then  becomes 


,,.4  I 


^  m=  I  > 3 , 5 . .  n=  I  ,3,5. .  p2m3r)+mn3  m+( p-s  1 


+  1  -  I  I 

m-(p-s )  m+(p+s)  m-(p- 


n-(q-i 


n+<q+tl 
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I 

n-(q++) 


(III- 11) 


where  i  corresponds  to  the  pair  of  integers  (p,  q),  and  k  corresponds  to 
the  pair  of  integers  (s,t).  For 


I  R 

=  J  J  4»^< ot , B )  da 


dB 


(III- 12) 


o  o 


then 


I  R 


i^(k,  i )  =  J  J  sin 


if n  q  sin  t  li  da  dB  <III_13) 


o  o 


p  it  a  sin  s  ira  sin  q  —  sin  t  —  da  dB 

R  R 


I 3<  k, i )  =  R/4  for  k  =  i 


=  0  for  k  ^  i 


(III- 14) 
(III- 15) 


For  the  special  case  of  the  square  duct  where  the  eigenfunctions 
used  are 


$j(a,6)  —  sin  p  ira  sin  q  ttB  +  sin  q  tto  sin  p  irB 


(III- 16) 


the  integrals  may  be  shown  to  be 

1  ( k,  i )  =  -  ir2(p2  +  q2)/2  k  =  i  and  p  i  q 


I  |  (k,  i )  =  -  it2 ( p2  +  q2)  k  =  i  and  p  =  q 


I  .  (k,  i )  =  0  k  l 


(III- 17) 
(III- 18) 

(III- 19) 
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where  i  corresponds  to  the  pair  of  integers  (p,  q),  and  k  corresponds 
to  the  pair  of  integers  (s,t). 


3<k,l)  =  1 

k  ■  l  and  p  =  q 

(III- 21) 

3(k,l)  =  1/2 

k  =  i  and  p  f  q 

(m-22) 

3<k,0  =  0 

k  *  i 

(III- 23) 
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APPENDIX  IV 

MISCELLANEOUS  INTEGRALS 


For  the  integrals  given  in  this  appendix  the  functions,  $i(a,  j3)  ,  are 
those  given  by  Eqs.  (26)  and  (27). 

The  integral  I(j£,m)  as  given  by  Eq.  (24), 


(4,m) 


R 

■-// 


$1 (a, 3) 


N 

L 

1=2 


C, 


(m) 


4> -  (a, 3) 


<frj(o,3)  + 


j=2 


do  d3 


(IV- 1) 


for  the  nonsquare  cases  becomes,  on  integration. 


<£,m)  =  - 

4 


+  £  Cj(m)  CjU) 


i =2 


(IV-2) 


and  for  the  square  cases  becomes 


IU,m)  = 


C* 


(IV-3) 


where 


C*  .=  Cj (m>  Cj <£)  for  p  =  q 


-I  r  (m)  r  (£)  ,  n  j  „ 

7  Li  Cj  for  p  t  q 
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The  integral  for  J\(i)  as  given  by  Eq.  (30), 

mN  1 


o  o  L 


j=2 


da  dB 


for  the  nonsquare  cases  becomes,  on  integration, 


Jx(£> 


r  n 

■S  '*  i 

L  j=2 


U) 


/pq 


and  for  the  square  cases  becomes 


V*>  • 


I  + 


I!  Cj^’/pq 


J-2 


The  integral  for  Jrj(£)  as  given  by  Eq.  (31), 

I  R 

v«>  -  //  L.  I  Ts(s,+)  sin  s  na  sin 

o  o  s= I ,3,5. .  t=l ,3,5. . 


A 

<pj  (a,&)  +  CjU)  <tij  (a, 6) 

j=2 


da  dB 


for  the  nonsquare  cases  becomes,  on  integration, 

N 


JnU>  =  4 


Ts  (1 


,d  +  £ 


f_(p,q)  C 


U> 


j=2 


j 


and  for  the  square  cases  becomes 


(IV-4) 


(IV- 5) 


(IV- 6) 


+  5JL 

R 


(IV -7) 


(IV-8) 
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JnU)  = 


N 


2  f g(  I ,  I )  +  2  fs(P'c’>  Cj 

j=2 


U> 


(IV -9) 


The  bulk  mean  temperature  which  is  defined  as 

I  R 


//  9(a,@,£)  u(a,B>  dad@ 


QJV  = 

m 


o  o 


I  R 


(IV- 10) 


J  J  u(a,8) 


da  dB 


o  o 


for  the  nonsquare  cases  becomes,  on  integration,  for  the  nonviscous 
component 


4RA0T 

0mxK)  =  — j-i-  y  a  (n) 
**  I A 

n=l 


R2+l 


N 


( n  > 


R2p3q+pq3 


-BnS/(IV-ll) 


4RPe 

»4  > 
n=l 


(n) 


R2+l 


IN 

I 

i -2 


(n)  Rz 


R2p3q+pq3 


Bnt. 

le  /  (IV- 12) 


and  for  the  square  cases  becomes,  for  the  nonviscous  component. 


_  I6A0t  V  (n)  I  .  1 

0mx(O  -  — —  2-  9I>.  4  +  2  L 


n=  I 


i  =2 


(n) 


p3q+pq3 


-Bny  (IV_  13) 


m 


and  for  the  viscous  component 

N 

mn  X.  iq 


v*  tn)  I  I  V 

L  '’in  4  +  2  Z- 


n=  I 


Ci 


(n)  1 


“  p  Vpq3  J 


“Bn* 


(IV- 14) 
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The  ratio  of  heat-transfer  rate  per  unit  length  at  the  duct  walls  to 
the  thermal  conductivity  may  be  written,  for  the  a  or  short  wall,  as 


VK  * 


(IV- 15) 


and  for  the  /3  or  long  wall  as 


VK  ■ 


/ 


3_0 

3a 


d6 

a=Q 


(IV- 16) 


and  becomes,  on  differentiation  and  integration,  for  the  nonsquare  cases 


N 

y  (n) 

Qa/  K  =  ai 

n=l 


N 


2*2  I 


i=2 


(n)  „ 

C:  9_ 

1  Rp 


e-Bn^ 


(IV- 17) 


VK 


■t 


(n) 


2R  +  2 


t 

i=2 


)  pR 

q 


e 


-Bne  (IV- 18) 


where  the  coefficients  ajW  become  a^fa)  anda^w  for  the  nonviscous 
and  viscous  components,  respectively,  and  for  the  square  cases  becomes 


vk-vk‘  Z 


(n) 


n=  I 


4  +  2 


N 

z 

i=2  » 


a+a 

q  p 


e"Bne  (IV- 19) 


with  the  above-mentioned  condition  on  the  coefficient  a . 

With  the  evaluation  of  the  heat-transfer  rates,  the  Nusselt  numbers 
may  now  be  written,  for  the  nonviscous  component,  as 


QaX  Q&X 

nr  +nc 


(R+l )emA 


(IV- 20) 
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and  for  the  viscous  component  as 


Nun 


_  h  L  | 
“  ~TC~ 
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